Time Step Size Research Articles

Nonsmooth Reduced Interface Models and Their Use in Co-Simulation of Mechanical Systems

Abstract In a co-simulation setup, the entire system is decomposed into a collection of individual subsystems that are interfaced together, with each subsystem being modeled and integrated separately according to its own requirements. To maintain the interconnectivity and consolidation of the primary system, these subsystems must communicate with each other through the interface and transfer certain information at the end points of a defined time interval termed macro time step. Inside the macro time step, the evolution of the interface variables has to be approximated as information about them will only be available again at the end of the step. In real-time simulations, the size of the macro time step and the accuracy of the approximated interface variables are critical factors; if the interface variables are approximated accurately, the size of the macro time step can be kept large enough to provide interactive rates without loss of accuracy and stability. This work focuses on systems where unilateral contact interactions are important and proposes reduced interface model concepts for such nonsmooth systems. The use of the proposed reduced interface model (RIM) is demonstrated in co-simulation to provide model-based approximation of the interface variables. The advantages of the proposed method are demonstrated through two representative case studies.

The Finite Element Method in Time for Multibody Dynamics

Abstract The generalized-α scheme has become the approach of choice for the time integration of the equations of motion of multibody systems. Despite its simplicity, this scheme presents drawbacks: the time-step size cannot be changed easily, making it difficult to implement time adaptivity, and the solution of periodic problems cannot be found easily. This paper explores an alternative approach based on the finite element method in time. The basic principles underpinning the approach are presented and both time-continuous and time-discontinuous approaches are investigated. Two types of Galerkin schemes will be presented here: the time-continuous and the time-discontinuous schemes. In the former, the displacement field is continuous across interelement boundaries, whereas discontinuities or “jumps” are allowed across interelement boundaries for the latter. Simple problems are treated to identify the best schemes. Families of schemes of various accuracies are presented. The first family, based on time-continuous elements, features schemes that do not present numerical dissipation. Asymptotic annihilation is achieved by the time-discontinuous elements that form the second family. The problem of kinematic constraints is treated within the framework of the finite element method in time. Special emphasis is devoted to the satisfaction of the kinematic constraints and their time derivative within a time element.

Dynamics of the Fractional-Order Lorenz System Based on Adomian Decomposition Method and Its DSP Implementation

Dynamics and digital circuit implementation of the fractional-order Lorenz system are investigated by employing Adomian decomposition method (ADM). Dynamics of the fractional-order Lorenz system with derivative order and parameter varying is analyzed by means of Lyapunov exponents (LEs), bifurcation diagram, chaos diagram and phase diagram. Results show that the fractional-order Lorenz system has rich dynamical behavior and it is a potential model for application. It is also found that the minimum order is affected by numerical algorithm and time step size. Finally, the fractional-order system is implemented on DSP digital circuit. Phase diagrams generated by the DSP are consistent with that generated by simulation.

Determination of Maneuvering Force Coefficients for a Destroyer Model with OpenFOAM

_ Ship maneuvering performance can be predicted using various methods, including physical model tests, rapid simulations using coefficient-based forces, and high-fidelity computational fluid dynamics. This article considers the application of computational fluid dynamics for the prediction of hull force coefficients to be used as inputs to rapid maneuvering simulations. The open-source software OpenFOAM was used to simulate forces on the destroyer model DTMB 5415 for steady drift, oscillatory pure sway motion, and oscillatory yaw motion. Recommendations are provided regarding best practices in a number of areas, including domain size, mesh refinement, time step size, and turbulent modeling. Predicted forces for steady drift angles up to 12 degrees are typically within 10%of experimental values. For oscillatory sway and yaw motions, predicted forces are typically within 25%of experimental values. Keywords maneuvering; OpenFOAM; best modeling practices

A GPU based accelerated solver for simulation of heat transfer during metal casting process

Abstract The metal casting process, which is one of the key drivers of the manufacturing industry, involves several physical phenomena occurring simultaneously like fluid flow, phase change, and heat transfer which affect the casting yield and quality. Casting process modeling involves numerical modeling of these phenomena on a computer. In recent decades, this has become an inevitable tool for foundry engineers to make defect-free castings. To expedite computational time Graphics Processing Units (GPUs) are being increasingly used in the numerical modeling of heat transfer and fluid flow. Initially, in this work a CPU based implicit solver code is developed for solving the 3D unsteady energy equation including phase change numerically using Finite Volume Method (FVM) which predicts the thermal profile during solidification in the metal casting process in a completely filled mold. To address the computational bottleneck, which is identified as the linear algebraic solver based on the Bi-Conjugate Gradient Stabilized (BiCGSTAB) method, a GPU-based code is developed using CUDA toolkit and was implemented on the GPU. The CPU and GPU based codes are then validated against a commercial casting simulation code FLOW-3D CAST® for a simple casting part and against in-house experimental results for gravity die casting of a simple geometry. Parallel performance is analyzed for grid sizes ranging from 10x10x10 to 210x210x210 and for three time-step sizes. The performance of the GPU code based on occupancy and throughput is also investigated. The GPU code exhibits a maximum speedup of 308x compared to the CPU code for a grid size of 210x210x210 and a time-step size of 2s.

Comprehensive investigation of the explicit, positivity preserving methods for the heat equation

In this paper-series, we investigate the performance of 12 explicit non-conventional algorithms. All of them have the convex combination property, thus they are unconditionally stable and preserve the positivity of the solution when they applied to the heat equation. In this part of the series, we construct several 2D systems to find how the errors depend on the time step size. Sweeps for other key parameters will be presented in the next part of the series.

Comprehensive investigation of the explicit, positivity preserving methods for the heat equation

In this paper-series, we investigate the performance of 12 explicit non-conventional algorithms in several 2D systems. All of them have the convex combination property, thus they are unconditionally stable and preserve the positivity of the solution when they are applied to the heat equation. In the first part of the series, we examined how the errors depend on the time step size and running times. Now we present additional numerical test results, where sweeps for parameters such as the stiffness and the wavelength of the initial function will be performed.

A highly parallel algorithm for simulating the elastodynamics of a patient-specific human heart with four chambers using a heterogeneous hyperelastic model

In this paper a highly parallel method is developed for simulating the elastodynamics of a four-chamber human heart with patient-specific geometry. The heterogeneous hyperelastic model is discretized by a finite element method in space and a fully implicit adaptive method in time, and the resulting nonlinear algebraic systems are solved by a scalable domain decomposition algorithm. The deformations of the cardiac muscles are quite complex due to the realistic geometry, the heterogeneous hyperelasticity of the cardiac tissue, and the myocardial fibers with active stresses. Moreover, the deformations in different chambers and at different phases of the cardiac cycle are very different. To simulate all the muscle movements including the atrial diastole, the atrial systole, the isovolumic contraction, the ventricular ejection, the isovolumic relaxation, and the ventricular filling, the temporal-spatial mesh needs to be sufficiently fine, but not too fine so that the overall computing time is manageable, we introduce a baseline mesh in space and a two-level time stepping strategy including a uniform baseline time step size to obtain the desired time accuracy and an adaptive time stepping method within a baseline time step to guarantee the convergence of the nonlinear solver. Through numerical experiments, we investigate the performance of the proposed method with respect to the material coefficients, the fiber orientations, as well as the mesh sizes and the time step sizes. For an unstructured tetrahedral mesh with more than 200 million degree of freedoms, the method scales well for up to 16,384 processor cores for all steps of an entire cardiac cycle.

One-Way Fluid–Structure Interaction Simulation of the Added Damping of T(0,1) Modes of Two Partially Overlapped Identical Plates in Still Water

Abstract Structures with a partially overlapped status in water can be seen in some engineering applications, and the fluid-structure coupling vibration behavior of two partially overlapped identical plates has been studied previously through the experimental method. In this study, the added damping of the T(0,1) modes of two submerged partially overlapped plates is numerically investigated through the one-way fluid-structure interaction (FSI) method. The relevant numerical settings, like the mesh size, calculated periods, time-step size, and vibration amplitude, were tested first. Then, the numerical results were compared with experimental results, and good agreements were found. Finally, numerical results were analyzed. The vibration status of two plates in the joint abutting area or the overlapped area has an important influence on the added mass variation. When the added mass is higher, the phase difference between modal force and vibration displacement is also greater, which is the main reason for the higher added damping. The relationship between the phase difference and the frequency in water can be approximately fitted to a straight line, which can probably be used to predict the added damping variations caused by fluid boundary changes of submerged structures.

Monolithic and local time-stepping decoupled algorithms for transport problems in fractured porous media

Abstract The objective of this paper is to develop efficient numerical algorithms for the linear advection-diffusion equation in fractured porous media. A reduced fracture model is considered where the fractures are treated as interfaces between subdomains and the interactions between the fractures and the surrounding porous medium are taken into account. The model is discretized by a backward Euler upwind-mixed hybrid finite element method in which the flux variable represents both the advective and diffusive fluxes. The existence, uniqueness, as well as optimal error estimates in both space and time for the fully discrete coupled problem are established. Moreover, to facilitate different time steps in the fracture-interface and the subdomains, global-in-time, nonoverlapping domain decomposition is utilized to derive two implicit iterative solvers for the discrete problem. The first method is based on the time-dependent Steklov–Poincaré operator, while the second one employs the optimized Schwarz waveform relaxation (OSWR) approach with Ventcel-Robin transmission conditions. A discrete space-time interface system is formulated for each method and is solved iteratively with possibly variable time step sizes. The convergence of the OSWR-based method with conforming time grids is also proved. Finally, numerical results in two dimensions are presented to verify the optimal order of convergence of the monolithic solver and to illustrate the performance of the two decoupled schemes with local time-stepping on problems of high Péclet numbers.

Development and benchmarking of space-time kinetics solver in REDAC

A space–time kinetics solver based on finite volume approach has been developed for REactor Dynamics Analysis Code (REDAC). It solves 3D multi-group neutron diffusion equation, coupled with the heat transport equations. Code offers flexibility in terms of reactor geometry, boundary conditions and discretization. Efficient numerical methods have been employed for solving the governing equations. The performance of the code has been demonstrated by solving diverse benchmark problems. It is shown that REDAC is able to predict important static and dynamic core parameters with good accuracy. The framework of REDAC easily lends itself for sensitivity studies with respect to spatial discretization and temporal integration schemes. From these studies, it is observed that spatial resolution could have significant influence on accuracy and the semi-implicit method performs better than the implicit method. Improved quasi-static method requires less computational time compared to the direct method, however, its accuracy depends on the time step size.

Enhanced physics‐informed neural networks for efficient modelling of hydrodynamics in river networks

AbstractThis article enhances the physics‐informed neural networks (PINNs) method to effectively model the hydrodynamics of real‐world river networks with irregular cross‐sections. First, we pre‐process hydraulic parameters to optimize training speed without compromising accuracy, achieving a 91.67% acceleration compared with traditional methods. To address the vanishing gradient problem, layer normalization is also incorporated into the architecture. We also introduce novel physical constraints—water level range and junction node equations—to ensure effective training convergence and enrich the model with additional physical insights. Two practical case studies using HEC‐RAS benchmarks demonstrate that our improved PINN method can predict river network hydrodynamics with less data and is less sensitive to time step size, allowing for longer computational time steps. Incorporating physical knowledge, our enhanced PINN methodology emerges as an efficient and promising avenue for modelling the complexities of hydrodynamic processes in natural river networks.

Efficient discrete element modeling of heat generation and transfer in granular flows: Validation and application

This work presents an efficient Discrete Element Method (DEM) framework for the simulation of the thermal behavior of granular media. The main focus is on long-lasting granular flows, involving heat transfer and generation due to mechanical energy dissipation. The proposed approach uses efficient strategies to reduce the computational cost of the analyses and, therefore, to enable its application to problems of practical relevance. For instance, the contact area is adjusted to compensate for the artificial material softening that is typically considered in DEM to increase the time step size. After extended validation, the methodology is applied to the simulation of different setups of an experimental rotating drum. The numerical simulations presented good agreement with the experimental results and allowed a detailed analysis of the mechanisms and patterns of heat generation, which could not be extrapolated from the experimental campaign.

Unconditional stability and error estimates of the modified characteristics FEMs for the Micropolar Navier–Stokes Equations

In this paper, the unconditional stability and optimal error estimate of the velocity, pressure and angular velocity for the modified characteristics FEMs of the unsteady Micropolar Naiver-Stokes Equations (MNSE) are presented. In this method, the nonlinear equation is linearized by the characteristic finite element method for dealing with the time derivative term and the convection term. Basing on the characteristic time-discrete system, the error function is split into a temporal error and a spatial error. With a rigorous analysis to the characteristic time-discrete system, we prove that the error between the numerical solution and the solution of the time-discrete system is τ-independent, where τ denotes the time stepsize. The stability results and optimal error estimates in L2 norm and H1 norm will be given. Finally, some numerical results will be provided to confirm our theoretical analysis.

Physics-Informed Neural Networks for Solving High-Index Differential-Algebraic Equation Systems Based on Radau Methods

As is well known, differential algebraic equations (DAEs), which are able to describe dynamic changes and underlying constraints, have been widely applied in engineering fields such as fluid dynamics, multi-body dynamics, mechanical systems, and control theory. In practical physical modeling within these domains, the systems often generate high-index DAEs. Classical implicit numerical methods typically result in varying order reduction of numerical accuracy when solving high-index systems. Recently, the physics-informed neural networks (PINNs) have gained attention for solving DAE systems. However, it faces challenges like the inability to directly solve high-index systems, lower predictive accuracy, and weaker generalization capabilities. In this paper, we propose a PINN computational framework, combined Radau IIA numerical method with an improved fully connected neural network structure, to directly solve high-index DAEs. Furthermore, we employ a domain decomposition strategy to enhance solution accuracy. We conduct numerical experiments with two classical high-index systems as illustrative examples, investigating how different orders and time-step sizes of the Radau IIA method affect the accuracy of neural network solutions. For different time-step sizes, the experimental results indicate that utilizing a 5th-order Radau IIA method in the PINN achieves a high level of system accuracy and stability. Specifically, the absolute errors for all differential variables remain as low as 10−6, and the absolute errors for algebraic variables are maintained at 10−5. Therefore, our method exhibits excellent computational accuracy and strong generalization capabilities, providing a feasible approach for the high-precision solution of larger-scale DAEs with higher indices or challenging high-dimensional partial differential algebraic equation systems.

A new linearized second-order energy-stable finite element scheme for the nonlinear Benjamin-Bona-Mahony-Burgers equation

In this paper, a novel linearized second-order energy-stable fully discrete scheme for the nonlinear Benjamin-Bona-Mahony-Burgers (BBMB) equation is introduced, along with a superconvergence analysis of the conforming finite element method (FEM). Through skillful decomposition of the nonlinear term u∇⋅u, a new linearized second-order fully discrete scheme is developed. This scheme requires fewer iterations and exhibits higher computational efficiency compared to the nonlinear approach. Furthermore, it maintains energy stability, enabling direct numerical solution boundedness in the H1-norm, which represents an improvement over the L∞-norm boundedness reported in previous literature. Secondly, relying on this boundedness and the high-precision results of the bilinear element, we derive the unconditional error estimates for superclose and global superconvergence without any restrictions on the time step size Δt and spatial mesh size h. Finally, numerical examples are presented to validate the theoretical analysis.

Simulation of nuclear fuel clad high-temperature ballooning under loss-of-coolant accident conditions considering anisotropic creep

Zirconium (Zr) alloys are standard clad material in commercial light water reactors. However, it is vulnerable to ballooning and burst under high-temperature conditions such as LOCA (Loss-Of-Coolant Accident), due to large internal/external pressure difference. Considering that clad ballooning is detrimental in terms of nuclear reactor safety especially for high-burnup fuel rod, high-fidelity numerical simulation of Zircaloy clad ballooning is desirable. This paper presented a finite-element based implementation for the simulation of Zircaloy clad ballooning. The developed simulation capability predicts progressive clad ballooning up to failure given empirical burst limits. The material models coupled together are phase transformation and high-temperature creep. Especially, anisotropic creep was considered for α phase Zircaloy. The implemented algorithm for stress update is implicit in time integration, which allows relatively large time step size to efficiently simulate the process that can last for hundreds of seconds. An analytic benchmark, the PUZRY ballooning and burst test, and the Hardy stress-rupture test were utilized to verify and validate the simulation capability. The anisotropic creep calculation conforms to the analytical solution. The comparison with PUZRY and Hardy test data shows that considering anisotropic creep in the α phase region is essential. 3D ballooning was also simulated by considering circumferential temperature difference for both conventional and novel annular clad designs. Interestingly, it is the first time to note that the risk of buckling failure is predicted for the double-cooled annular nuclear fuel clad during the accident condition, only after considering the anisotropic creep algorithm.

Learning Reduced Fluid Dynamics

Predicting the state evolution of ultra high-dimensional, time-reversible fluid dynamic systems is a crucial but computationally expensive task. Existing physics-informed neural networks either incur high inference cost or cannot preserve the time-reversible nature of the underlying dynamics system. We propose a model-based approach to identify low-dimensional, time reversible, nonlinear fluid dynamic systems. Our method utilizes the symplectic structure of reduced Eulerian fluid and use stochastic Riemann optimization to obtain a low-dimensional bases that minimize the expected trajectory-wise dimension-reduction error over a given distribution of initial conditions. We show that such minimization is well-defined since the reduced trajectories are differentiable with respect to the subspace bases over the entire Grassmannian manifold, under proper choices of timestep sizes and numerical integrators. Finally, we propose a loss function measuring the trajectory-wise discrepancy between the original and reduced models. By tensor precomputation, we show that gradient information of such loss function can be evaluated efficiently over a long trajectory without time-integrating the high-dimensional dynamic system. Through evaluations on a row of simulation benchmarks, we show that our method reduces the discrepancy by 50-90 percent over conventional reduced models and we outperform PINNs by exactly preserving the time reversibility.

Effect of Different Time Step Sizes on Pedestrian Evacuation Time under Emergencies Such as Fires Using an Extended Cellular Automata Model

The cellular automata (CA) model has been a meaningful way to study pedestrian evacuation during emergencies, such as fires, for many years. Although the time step used in the CA model is one of the most essential elements, there is a lack of research on its impact on evacuation time. In this paper, we set different time step sizes in an extended cellular automaton model and discuss the effect of time step size on the overall evacuation time under different emergency types and levels. For a fixed step time mode, the larger the time step, the longer the evacuation time. In each time step size, the evacuation time gradually increases with the increase of emergency level, and there is a sharp increase when the time for pedestrians to move one step is exactly an integer multiple of the time step. When there is no friction between pedestrians, the evacuation time at each time step first decreases slightly with the increase of emergency level and then remains unchanged; the larger the time step, when the evacuation time remains unchanged, the lower the emergency level and the greater the evacuation time. For the variable time step model, when the friction between pedestrians approaches infinity, the total evacuation time does not change with the emergency level; when the friction between pedestrians is reduced, the total evacuation time slightly decreases with the increase of the emergency level. The less friction there is, the more significant the reduction. The results of previous actual experiments are also reflected in the simulation at a lower emergency level. The result shows that the time step size significantly impacts the evacuation simulation results of the CA model, and researchers should choose carefully to obtain more realistic simulation results.

A time two-grid algorithm for two-dimensional nonlinear time-fractional partial integro-differential equations

In this paper, a temporal second order two-grid difference scheme is proposed for the two-dimensional nonlinear time-fractional partial integro-differential equations with a weakly singular kernel. The first-order backward difference and L1 formula are used in the temporal direction to estimate the first level of time, the L2−1σ formula and L1-type formula are used in the temporal direction for later time steps, and the central difference formula is used in the spatial directions. To improve the computational efficiency of nonlinear system, an efficient time two-grid algorithm is proposed. This algorithm firstly solves a nonlinear system on the coarse grid, and then the Lagrangian linear interpolation is applied on the coarse grid to estimate the function values on the fine grid. The stability and convergence of the two-grid difference scheme are analyzed by the energy method. The convergence order of the two-grid difference scheme is O(τF2+τC4+hx2+hy2), where τF and τC are the time step sizes of fine grid and coarse grid respectively, while hx and hy are the space step sizes. Numerical experiments show that the accuracy of the theoretical analysis and the efficiency of the two-grid algorithm.